function [solzarray,xarray]= provideFamilySolution_matrixnormPOLYFMIN(F,w,h,numparams,numsols,Kguess,type)
if(numparams~=2)
    error('only works with 2');
end

if(exist('type','var')~=1)
    type='foc1';
end

if(size(F,1)==0)
    solzarray=[];
    xarray=[];
    
    return ;
end

focalscale=3;

if(strcmpi(type, 'foc1')==1)
    [solzarray,xarray,xfarray]=findFocal11(F,w,h,numsols,Kguess,focalscale);
elseif(strcmpi(type, 'foc2')==1)
    [solzarray,xarray,xfarray]=findFocal12(F,w,h,numsols,Kguess,focalscale);
elseif(strcmpi(type, 'xc')==1)
    [solzarray,xarray,xfarray]=findOCxc(F,w,h,numsols,Kguess);
elseif(strcmpi(type, 'yc')==1)
    [solzarray,xarray,xfarray]=findOCyc(F,w,h,numsols,Kguess);
else
    error('no such option');
end

if( size(xfarray,1)==size(xarray,1))
    xarray=xfarray;
elseif(size(xfarray,1)>10 && size(xfarray,1)~=size(xarray,1))
    solzarray = interp1q(xfarray,solzarray,xarray);
    solzarray(isnan(solzarray))=0;
elseif(size(xfarray,1)~=size(xarray,1))
    solzarray=[];
    xarray=[];
end


end




function [solzarray,xarray,xfarray]=findFocal11(F,w,h,numsols,Kguess,focalscale)


F=F/norm(F,'fro');

focal2start=10;
focal2end=focalscale*w;

if(numsols==1)
    xarray=[Kguess(2,2)];
else
    step=((focal2end-focal2start)/(numsols-1));
    xarray=focal2start:step:focal2end;
    
end

solzarray=[];
xfarray=[];

f1=F(1,1);
f2=F(1,2);
f3=F(1,3);

f4=F(2,1);
f5=F(2,2);
f6=F(2,3);

f7=F(3,1);
f8=F(3,2);
f9=F(3,3);


xc=Kguess(1,3);
yc=Kguess(2,3);



minans=1;
maxans=focal2end;

exp1=f1^4;

exp2_a=(4*f1^2*f2^2 - 2*f1^2*(f2^2 + f4^2) + 4*f1^2*f4^2);
exp2_b= 4*f1^2*(f1^2*xc^2 + 2*f1*f2*xc*yc + 2*f1*f3*xc + f2^2*yc^2 + 2*f2*f3*yc + f3^2) + 4*f1^2*(f7 + f1*xc + f4*yc)^2 - 2*f1^2*((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2);

exp3_a=(2*f2^4 - (f2^2 + f4^2)^2 + 2*f4^4 - 2*f1^2*f5^2 + 8*f1*f2*f4*f5);
exp3_b= (4*f2^2*(f1^2*xc^2 + 2*f1*f2*xc*yc + 2*f1*f3*xc + f2^2*yc^2 + 2*f2*f3*yc + f3^2) + 4*f4^2*(f7 + f1*xc + f4*yc)^2 - 2*f1^2*((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2) - 2*((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2)*(f2^2 + f4^2) + 8*f1*f4*(f3*f6 + f1*f4*xc^2 + f2*f5*yc^2 + f1*f6*xc + f3*f4*xc + f2*f6*yc + f3*f5*yc + f1*f5*xc*yc + f2*f4*xc*yc) + 8*f1*f2*(f7 + f1*xc + f4*yc)*(f8 + f2*xc + f5*yc));
exp3_c=- ((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2)^2 + 2*(f7 + f1*xc + f4*yc)^4 + 2*(f1^2*xc^2 + 2*f1*f2*xc*yc + 2*f1*f3*xc + f2^2*yc^2 + 2*f2*f3*yc + f3^2)^2 - 2*f1^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 8*f1*(f3 + f1*xc + f2*yc)*(f7 + f1*xc + f4*yc)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc);

exp4_a= (4*f2^2*f5^2 - f5^2*(2*f2^2 + 2*f4^2) + 4*f4^2*f5^2);
exp4_b= (- f5^2*(2*(f3 + f1*xc + f2*yc)^2 + 2*(f7 + f1*xc + f4*yc)^2) + 4*f4^2*(f4^2*xc^2 + 2*f4*f5*xc*yc + 2*f4*f6*xc + f5^2*yc^2 + 2*f5*f6*yc + f6^2) - ((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2)*(2*f2^2 + 2*f4^2) + 4*f2^2*(f8 + f2*xc + f5*yc)^2 + 8*f2*f5*(f3*f6 + f1*f4*xc^2 + f2*f5*yc^2 + f1*f6*xc + f3*f4*xc + f2*f6*yc + f3*f5*yc + f1*f5*xc*yc + f2*f4*xc*yc) + 8*f4*f5*(f7 + f1*xc + f4*yc)*(f8 + f2*xc + f5*yc)) ;
exp4_c=(4*(f7 + f1*xc + f4*yc)^2*(f8 + f2*xc + f5*yc)^2 - (2*(f3 + f1*xc + f2*yc)^2 + 2*(f7 + f1*xc + f4*yc)^2)*((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2) + 4*(f3*f6 + f1*f4*xc^2 + f2*f5*yc^2 + f1*f6*xc + f3*f4*xc + f2*f6*yc + f3*f5*yc + f1*f5*xc*yc + f2*f4*xc*yc)^2 - (2*f2^2 + 2*f4^2)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 8*f2*(f3 + f1*xc + f2*yc)*(f8 + f2*xc + f5*yc)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc) + 8*f4*(f7 + f1*xc + f4*yc)*(f6 + f4*xc + f5*yc)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)) ;
exp4_d=4*(f3 + f1*xc + f2*yc)^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 - (2*(f3 + f1*xc + f2*yc)^2 + 2*(f7 + f1*xc + f4*yc)^2)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 4*(f7 + f1*xc + f4*yc)^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 ;

exp6_a=f5^4;
exp6_b=(4*f5^2*(f4^2*xc^2 + 2*f4*f5*xc*yc + 2*f4*f6*xc + f5^2*yc^2 + 2*f5*f6*yc + f6^2) + 4*f5^2*(f8 + f2*xc + f5*yc)^2 - 2*f5^2*((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2)) ;
exp6_c= (- ((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2)^2 + 2*(f8 + f2*xc + f5*yc)^4 + 2*(f4^2*xc^2 + 2*f4*f5*xc*yc + 2*f4*f6*xc + f5^2*yc^2 + 2*f5*f6*yc + f6^2)^2 - 2*f5^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 8*f5*(f6 + f4*xc + f5*yc)*(f8 + f2*xc + f5*yc)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)) ;
exp6_d= (4*(f6 + f4*xc + f5*yc)^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 4*(f8 + f2*xc + f5*yc)^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 - 2*((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2) ;
exp6_e=  (f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^4;

for i=1:numsols
    
    foc2=xarray(i);
    
    foc2p2=foc2^2;
    foc2p4=foc2p2^2;
    foc2p6=foc2p2^3;
    foc2p8=foc2p4^2;
    
    p=[ exp1
        0
        (exp2_a*foc2p2)+exp2_b
        0
        (exp3_a*foc2p4)+(exp3_b*foc2p2)+(exp3_c)
        0
        (exp4_a*foc2p6)+(exp4_b*foc2p4)+(exp4_c*foc2p2)+exp4_d
        0
        (exp6_a*foc2p8)+(exp6_b*foc2p6)+(exp6_c*foc2p4)+(exp6_d*foc2p2)+exp6_e ];
    
   f = @(x)((p(1)*x.^8)+(p(2)*x.^7)+(p(3)*x.^6)+(p(4)*x.^5)+(p(5)*x.^4)+(p(6)*x.^3)+(p(7)*x.^2)+(p(8)*x.^1)+(p(9)*x.^0));
    answ=processRoots_focalonly(fminbnd(f, focal2start, focal2end));
  
    if(answ>minans && answ<maxans)
        solzarray=[solzarray ;answ];
        xfarray= [xfarray ; xarray(i)];
        
    end
    
end
xarray=xarray';



end


function [solzarray,xarray,xfarray]=findFocal12(F,w,h,numsols,Kguess,focalscale)


F=F/norm(F,'fro');

focal1start=10;
focal1end=focalscale*w;
if(numsols==1)
    xarray=[Kguess(2,2)];
else
    step=((focal1end-focal1start)/(numsols-1));
    xarray=focal1start:step:focal1end;
end

solzarray=[];
xfarray=[];



f1=F(1,1);
f2=F(1,2);
f3=F(1,3);

f4=F(2,1);
f5=F(2,2);
f6=F(2,3);

f7=F(3,1);
f8=F(3,2);
f9=F(3,3);


xc=Kguess(1,3);
yc=Kguess(2,3);

maxans=focal1end;
minans=1;

exp1=f5^4 ;

exp2_a= (4*f2^2*f5^2 - 2*f5^2*(f2^2 + f4^2) + 4*f4^2*f5^2);
exp2_b=  4*f5^2*(f4^2*xc^2 + 2*f4*f5*xc*yc + 2*f4*f6*xc + f5^2*yc^2 + 2*f5*f6*yc + f6^2) + 4*f5^2*(f8 + f2*xc + f5*yc)^2 - 2*f5^2*((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2) ;

exp3_a= (2*f2^4 - (f2^2 + f4^2)^2 + 2*f4^4 - 2*f1^2*f5^2 + 8*f1*f2*f4*f5);
exp3_b= (4*f4^2*(f4^2*xc^2 + 2*f4*f5*xc*yc + 2*f4*f6*xc + f5^2*yc^2 + 2*f5*f6*yc + f6^2) + 4*f2^2*(f8 + f2*xc + f5*yc)^2 - 2*f5^2*((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2) - 2*((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2)*(f2^2 + f4^2) + 8*f2*f5*(f3*f6 + f1*f4*xc^2 + f2*f5*yc^2 + f1*f6*xc + f3*f4*xc + f2*f6*yc + f3*f5*yc + f1*f5*xc*yc + f2*f4*xc*yc) + 8*f4*f5*(f7 + f1*xc + f4*yc)*(f8 + f2*xc + f5*yc)) ;
exp3_c= - ((f6 + f4*xc + f5*yc)^2 + (f8 + f2*xc + f5*yc)^2)^2 + 2*(f8 + f2*xc + f5*yc)^4 + 2*(f4^2*xc^2 + 2*f4*f5*xc*yc + 2*f4*f6*xc + f5^2*yc^2 + 2*f5*f6*yc + f6^2)^2 - 2*f5^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 8*f5*(f6 + f4*xc + f5*yc)*(f8 + f2*xc + f5*yc)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc);

exp4_a= (4*f1^2*f2^2 - f1^2*(2*f2^2 + 2*f4^2) + 4*f1^2*f4^2) ;
exp4_b=  (- f1^2*(2*(f6 + f4*xc + f5*yc)^2 + 2*(f8 + f2*xc + f5*yc)^2) + 4*f2^2*(f1^2*xc^2 + 2*f1*f2*xc*yc + 2*f1*f3*xc + f2^2*yc^2 + 2*f2*f3*yc + f3^2) - ((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2)*(2*f2^2 + 2*f4^2) + 4*f4^2*(f7 + f1*xc + f4*yc)^2 + 8*f1*f4*(f3*f6 + f1*f4*xc^2 + f2*f5*yc^2 + f1*f6*xc + f3*f4*xc + f2*f6*yc + f3*f5*yc + f1*f5*xc*yc + f2*f4*xc*yc) + 8*f1*f2*(f7 + f1*xc + f4*yc)*(f8 + f2*xc + f5*yc)) ;
exp4_c=  (4*(f7 + f1*xc + f4*yc)^2*(f8 + f2*xc + f5*yc)^2 - (2*(f6 + f4*xc + f5*yc)^2 + 2*(f8 + f2*xc + f5*yc)^2)*((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2) + 4*(f3*f6 + f1*f4*xc^2 + f2*f5*yc^2 + f1*f6*xc + f3*f4*xc + f2*f6*yc + f3*f5*yc + f1*f5*xc*yc + f2*f4*xc*yc)^2 - (2*f2^2 + 2*f4^2)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 8*f2*(f3 + f1*xc + f2*yc)*(f8 + f2*xc + f5*yc)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc) + 8*f4*(f7 + f1*xc + f4*yc)*(f6 + f4*xc + f5*yc)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)) ;
exp4_d=  4*(f6 + f4*xc + f5*yc)^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 - (2*(f6 + f4*xc + f5*yc)^2 + 2*(f8 + f2*xc + f5*yc)^2)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 4*(f8 + f2*xc + f5*yc)^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 ;

exp6_a=f1^4 ;
exp6_b=  (4*f1^2*(f1^2*xc^2 + 2*f1*f2*xc*yc + 2*f1*f3*xc + f2^2*yc^2 + 2*f2*f3*yc + f3^2) + 4*f1^2*(f7 + f1*xc + f4*yc)^2 - 2*f1^2*((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2)) ;
exp6_c=  (- ((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2)^2 + 2*(f7 + f1*xc + f4*yc)^4 + 2*(f1^2*xc^2 + 2*f1*f2*xc*yc + 2*f1*f3*xc + f2^2*yc^2 + 2*f2*f3*yc + f3^2)^2 - 2*f1^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 8*f1*(f3 + f1*xc + f2*yc)*(f7 + f1*xc + f4*yc)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc));
exp6_d=   (4*(f3 + f1*xc + f2*yc)^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 + 4*(f7 + f1*xc + f4*yc)^2*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2 - 2*((f3 + f1*xc + f2*yc)^2 + (f7 + f1*xc + f4*yc)^2)*(f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^2);
exp6_e=  (f9 + f3*xc + f7*xc + f6*yc + f8*yc + f1*xc^2 + f5*yc^2 + f2*xc*yc + f4*xc*yc)^4  ;


for i=1:numsols
    
    foc1=xarray(i);
    
    foc1p2=foc1^2;
    foc1p4=foc1p2^2;
    foc1p6=foc1p2^3;
    foc1p8=foc1p4^2;
    
    
    p=[ exp1
        0
        (exp2_a*foc1p2)+exp2_b
        0
        (exp3_a*foc1p4)+(exp3_b*foc1p2)+(exp3_c)
        0
        (exp4_a*foc1p6)+(exp4_b*foc1p4)+(exp4_c*foc1p2)+exp4_d
        0
        (exp6_a*foc1p8)+(exp6_b*foc1p6)+(exp6_c*foc1p4)+(exp6_d*foc1p2)+exp6_e];
    
     f = @(x)((p(1)*x.^8)+(p(2)*x.^7)+(p(3)*x.^6)+(p(4)*x.^5)+(p(5)*x.^4)+(p(6)*x.^3)+(p(7)*x.^2)+(p(8)*x.^1)+(p(9)*x.^0));
    answ=processRoots_focalonly(fminbnd(f, focal1start, focal1end));
    
    if(answ>minans && answ<maxans)
        solzarray=[solzarray ;answ];
        xfarray= [xfarray ; xarray(i)];
        
    end
    
end
xarray=xarray';


end


function [solzarray,xarray,xfarray]=findOCxc(F,w,h,numsols,Kguess)

F=F/norm(F,'fro');



xcstart=10;
xcend=w;
if(numsols==1)
    xarray=[Kguess(2,2)];
else
    step=((xcend-xcstart)/(numsols-1));
    xarray=xcstart:step:xcend;
    
end

solzarray=[];
xfarray=[];



maxans=xcend;
minans=1;

f1=F(1,1);
f2=F(1,2);
f3=F(1,3);

f4=F(2,1);
f5=F(2,2);
f6=F(2,3);

f7=F(3,1);
f8=F(3,2);
f9=F(3,3);

foc1=Kguess(1,1);
foc2=Kguess(2,2);



exp1= f1^4;

exp2_a=(4*f1^3*(f2 + f4)) ;
exp2_b= 4*f1^3*(f3 + f7);

exp3_a= (2*f1^2*(2*f1*f5 + (f2 + f4)^2) + 4*f1^2*(f2 + f4)^2);
exp3_b=  (2*f1^2*(2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8)) + 8*f1^2*(f2 + f4)*(f3 + f7));
exp3_c=  4*f1^2*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 4*f1^2*(f3 + f7)^2 + 4*f1^4*foc1^2 - 2*f1^2*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 4*f1^2*f4^2*foc2^2;


exp4_a= (4*f1*(2*f1*f5 + (f2 + f4)^2)*(f2 + f4) + 2*f1^2*f5*(2*f2 + 2*f4));
exp4_b=  (4*f1*((2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(f2 + f4) + (2*f1*f5 + (f2 + f4)^2)*(f3 + f7)) + 2*f1^2*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7)));
exp4_c= (- 4*f1*((f2 + f4)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(f3 + f7)) + 8*f1*((f2 + f4)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(f3 + f7)) - 2*f1^2*((f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2) + 4*f1^2*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4)) + 8*f1^2*foc1^2*(f1*f2 + f1*(f2 + f4)) + 8*f1*f4*foc2^2*(f1*f5 + f4*(f2 + f4))) ;
exp4_d=  4*f1^2*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + f9*(2*f3 + 2*f7)) - 2*f1^2*(f9*(2*f3 + 2*f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 8*f1^2*foc1^2*(f1*f3 + f1*(f3 + f7)) - 4*f1*(f3 + f7)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 8*f1*(f3 + f7)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 8*f1*f4*foc2^2*(f1*f6 + f4*(f3 + f7));

exp5_a= ((2*f1*f5 + (f2 + f4)^2)^2 + 2*f1^2*f5^2 + 4*f1*f5*(f2 + f4)*(2*f2 + 2*f4));
exp5_b=  (4*f1*((f2 + f4)*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7)) + f5*(f3 + f7)*(2*f2 + 2*f4)) + 2*(2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(2*f1*f5 + (f2 + f4)^2) + 4*f1^2*f5*(f6 + f8));
exp5_c= (4*foc1^2*(f1*f2 + f1*(f2 + f4))^2 + 4*foc2^2*(f1*f5 + f4*(f2 + f4))^2 - 4*f1*((f3 + f7)*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7)) + (f2 + f4)*((f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2)) + 4*f1^2*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 8*f1*((f3 + f7)*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7)) + (f2 + f4)*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4))) - 2*(2*f1*f5 + (f2 + f4)^2)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 4*(2*f1*f5 + (f2 + f4)^2)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))^2 - 2*f1^2*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 8*f1^2*foc1^2*(f1*f5 + f2*(f2 + f4)) + 8*f1*f4*foc2^2*(f4*f5 + f5*(f2 + f4))) ;
exp5_d=  (- 2*f1^2*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) - 2*(2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 4*(2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 4*f1^2*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8)) - 4*f1*((f2 + f4)*(f9*(2*f3 + 2*f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + (f3 + f7)*((f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2)) + 8*f1*((f2 + f4)*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + f9*(2*f3 + 2*f7)) + (f3 + f7)*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4))) + 8*foc1^2*(f1*f2 + f1*(f2 + f4))*(f1*f3 + f1*(f3 + f7)) + 8*foc2^2*(f1*f5 + f4*(f2 + f4))*(f1*f6 + f4*(f3 + f7)) + 8*f1^2*foc1^2*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)) + 8*f1*f4*foc2^2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)));
exp5_e=  (4*foc1^2*(f1*f3 + f1*(f3 + f7))^2 + 4*foc2^2*(f1*f6 + f4*(f3 + f7))^2 + 4*f1^2*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 2*f1^2*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2)^2 + 2*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2)^2 + 2*f1^4*foc1^4 + 2*f4^4*foc2^4 + 8*f1^2*foc1*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) - 4*f1*(f3 + f7)*(f9*(2*f3 + 2*f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 8*f1*(f3 + f7)*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + f9*(2*f3 + 2*f7)) + 4*f1^2*f4^2*foc1^2*foc2^2 + 8*f1*f4*foc2*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2));

exp6_a= (4*f1*f5^2*(f2 + f4) + 2*f5*(2*f1*f5 + (f2 + f4)^2)*(2*f2 + 2*f4));
exp6_b=  (4*f1*((f3 + f7)*f5^2 + 2*(f2 + f4)*(f6 + f8)*f5) - (2*f1*f5 + (f2 + f4)^2)*(2*(f6 + f8)*(2*f2 + 2*f4) + 2*f5*(2*f3 + 2*f7)) + (2*f1*f5 + (f2 + f4)^2)*(4*(f6 + f8)*(2*f2 + 2*f4) + 4*f5*(2*f3 + 2*f7)) + 2*f5*(2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(2*f2 + 2*f4));
exp6_c=  (- 4*f1*((f2 + f4)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 2*f5*(f3 + f7)*(f6 + f8)) + 8*f1*((f2 + f4)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 2*f5*(f3 + f7)*(f6 + f8)) - (2*f1*f5 + (f2 + f4)^2)*(2*(f6 + f8)*(2*f3 + 2*f7) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + (2*f1*f5 + (f2 + f4)^2)*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f6 + f8)*(2*f3 + 2*f7) + 4*f9*(2*f2 + 2*f4)) - (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(2*(f6 + f8)*(2*f2 + 2*f4) + 2*f5*(2*f3 + 2*f7)) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(4*(f6 + f8)*(2*f2 + 2*f4) + 4*f5*(2*f3 + 2*f7)) + 8*foc1^2*(f1*f2 + f1*(f2 + f4))*(f1*f5 + f2*(f2 + f4)) + 8*foc2^2*(f1*f5 + f4*(f2 + f4))*(f4*f5 + f5*(f2 + f4)) - 2*f5*(2*f2 + 2*f4)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 4*f5*(2*f2 + 2*f4)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 8*f1^2*f2*f5*foc1^2 + 8*f1*f4*f5^2*foc2^2);
exp6_d=  (- 4*f1*((f3 + f7)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + (f2 + f4)*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)) - (2*f1*f5 + (f2 + f4)^2)*(2*f9*(2*f3 + 2*f7) + 4*f1*f3*foc1^2 + 4*f1*f7*foc1^2 + 4*f2*f8*foc2^2 + 4*f4*f6*foc2^2) + 8*foc1*(foc1*(f1*f2 + f1*(f2 + f4))*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)) + foc1*(f1*f5 + f2*(f2 + f4))*(f1*f3 + f1*(f3 + f7))) + 8*foc2*(foc2*(f1*f5 + f4*(f2 + f4))*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)) + foc2*(f4*f5 + f5*(f2 + f4))*(f1*f6 + f4*(f3 + f7))) + 8*f1*((f3 + f7)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + (f2 + f4)*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8))) + (2*f1*f5 + (f2 + f4)^2)*(8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 4*f9*(2*f3 + 2*f7)) - (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(2*(f6 + f8)*(2*f3 + 2*f7) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) - (2*(f6 + f8)*(2*f2 + 2*f4) + 2*f5*(2*f3 + 2*f7))*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f6 + f8)*(2*f3 + 2*f7) + 4*f9*(2*f2 + 2*f4)) + (4*(f6 + f8)*(2*f2 + 2*f4) + 4*f5*(2*f3 + 2*f7))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 8*f1^2*foc1^2*(f3*f5 + f2*(f6 + f8)) + 8*f1*f4*foc2^2*(f5*f6 + f5*(f6 + f8)));
exp6_e= (- (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(2*f9*(2*f3 + 2*f7) + 4*f1*f3*foc1^2 + 4*f1*f7*foc1^2 + 4*f2*f8*foc2^2 + 4*f4*f6*foc2^2) + 8*foc1*((f1*f2 + f1*(f2 + f4))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + foc1*(f1*f3 + f1*(f3 + f7))*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 4*f9*(2*f3 + 2*f7)) - 4*f1*((f2 + f4)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (f3 + f7)*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)) + 8*foc2*((f1*f5 + f4*(f2 + f4))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + foc2*(f1*f6 + f4*(f3 + f7))*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))) + 8*f1*((f2 + f4)*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (f3 + f7)*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8))) - (2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2)*(2*(f6 + f8)*(2*f3 + 2*f7) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + (8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f6 + f8)*(2*f3 + 2*f7) + 4*f9*(2*f2 + 2*f4))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 8*f1^3*f2*foc1^4 + 8*f4^3*f5*foc2^4 + 8*f1^2*foc1*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + 8*f1*f4*foc2*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + 8*f1*f4*foc1^2*foc2^2*(f1*f5 + f2*f4));
exp6_f=  - ((2*f9*(2*f3 + 2*f7) + 4*f1*f3*foc1^2 + 4*f1*f7*foc1^2 + 4*f2*f8*foc2^2 + 4*f4*f6*foc2^2)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) - (8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 4*f9*(2*f3 + 2*f7))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) - 8*foc1*(f1*f3 + f1*(f3 + f7))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) - 8*foc2*(f1*f6 + f4*(f3 + f7))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) - 8*f1^3*f3*foc1^4 - 8*f4^3*f6*foc2^4 - 8*f1^2*foc1*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) - 8*f1*(f3 + f7)*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 4*f1*(f3 + f7)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 8*f1*f4*foc2*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) - 8*f1*f4*foc1^2*foc2^2*(f1*f6 + f3*f4));

exp7_a= (4*f5^2*(2*f1*f5 + (f2 + f4)^2) + f5^2*(2*f2 + 2*f4)^2 - f5^2*(4*f1*f5 + 2*(f2 + f4)^2));
exp7_b= (4*f5^2*(2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8)) - f5^2*(4*(f2 + f4)*(f3 + f7) + 4*f1*(f6 + f8)) + 8*f5*(2*f1*f5 + (f2 + f4)^2)*(f6 + f8) + 2*f5*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7))*(2*f2 + 2*f4) - 2*f5*(f6 + f8)*(4*f1*f5 + 2*(f2 + f4)^2));
exp7_c= (4*foc1^2*(f1*f5 + f2*(f2 + f4))^2 + 4*foc2^2*(f4*f5 + f5*(f2 + f4))^2 + 4*f5^2*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) - f5^2*(4*f1*f9 + 2*(f3 + f7)^2 + 4*f1^2*foc1^2 + 2*f2^2*foc2^2 + 2*f4^2*foc2^2) - (4*f1*f5 + 2*(f2 + f4)^2)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + (2*f1*f5 + (f2 + f4)^2)*(8*f5*f9 + 4*(f6 + f8)^2 + 4*f4^2*foc1^2 + 4*f5^2*foc2^2) + ((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7))^2 + 4*f5*(2*f2 + 2*f4)*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4)) + 8*f5^2*foc2^2*(f1*f5 + f4*(f2 + f4)) + 8*f5*(2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(f6 + f8) - 2*f5*(4*(f2 + f4)*(f3 + f7) + 4*f1*(f6 + f8))*(f6 + f8) - 2*f5*(2*f2 + 2*f4)*((f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2) + 8*f2*f5*foc1^2*(f1*f2 + f1*(f2 + f4)));
exp7_d= (- 2*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7))*((f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2) - (4*(f2 + f4)*(f3 + f7) + 4*f1*(f6 + f8))*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 4*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7))*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4)) + (2*f1*f5 + (f2 + f4)^2)*(8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 8*f9*(f6 + f8)) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(8*f5*f9 + 4*(f6 + f8)^2 + 4*f4^2*foc1^2 + 4*f5^2*foc2^2) + 8*foc2*(f5^2*foc2*(f1*f6 + f4*(f3 + f7)) + foc2*(f1*f5 + f4*(f2 + f4))*(f5*f6 + f5*(f6 + f8))) - (4*f1*f5 + 2*(f2 + f4)^2)*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 8*foc1*(foc1*(f1*f2 + f1*(f2 + f4))*(f3*f5 + f2*(f6 + f8)) + f2*f5*foc1*(f1*f3 + f1*(f3 + f7))) - 2*f5*(f6 + f8)*(4*f1*f9 + 2*(f3 + f7)^2 + 4*f1^2*foc1^2 + 2*f2^2*foc2^2 + 2*f4^2*foc2^2) - 2*f5*(2*f2 + 2*f4)*(f9*(2*f3 + 2*f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 4*f5*(2*f2 + 2*f4)*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + f9*(2*f3 + 2*f7)) + 8*foc1^2*(f1*f5 + f2*(f2 + f4))*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)) + 8*foc2^2*(f4*f5 + f5*(f2 + f4))*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)) + 8*f5*(f6 + f8)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2));
exp7_e=  (- (4*f1*f9 + 2*(f3 + f7)^2 + 4*f1^2*foc1^2 + 2*f2^2*foc2^2 + 2*f4^2*foc2^2)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) - (4*(f2 + f4)*(f3 + f7) + 4*f1*(f6 + f8))*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + (8*f5*f9 + 4*(f6 + f8)^2 + 4*f4^2*foc1^2 + 4*f5^2*foc2^2)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 8*f9*(f6 + f8)) + 8*foc2*((f1*f5 + f4*(f2 + f4))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + foc2*(f1*f6 + f4*(f3 + f7))*(f5*f6 + f5*(f6 + f8))) - ((f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2)^2 + 8*foc1*((f1*f2 + f1*(f2 + f4))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + foc1*(f1*f3 + f1*(f3 + f7))*(f3*f5 + f2*(f6 + f8))) + 4*foc1^2*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))^2 + 4*foc2^2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))^2 - 2*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7))*(f9*(2*f3 + 2*f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 2*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4))^2 + (2*f1*f5 + (f2 + f4)^2)*(4*f7^2*foc1^2 + 4*f8^2*foc2^2 + 4*f9^2) + 4*((f6 + f8)*(2*f2 + 2*f4) + f5*(2*f3 + 2*f7))*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + f9*(2*f3 + 2*f7)) - (4*f1*f5 + 2*(f2 + f4)^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 8*foc1*(f1*f5 + f2*(f2 + f4))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + 8*foc2*(f4*f5 + f5*(f2 + f4))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + 4*foc1^2*foc2^2*((f1*f5 + f2*f4)^2 + 2*f1*f2*f4*f5) + 12*f1^2*f2^2*foc1^4 + 12*f4^2*f5^2*foc2^4);
exp7_f=  (8*foc1*((f1*f2 + f1*(f2 + f4))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + (f1*f3 + f1*(f3 + f7))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1)) - (4*(f2 + f4)*(f3 + f7) + 4*f1*(f6 + f8))*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 8*f9*(f6 + f8))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) - 2*(f9*(2*f3 + 2*f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)*((f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2) - (2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)*(4*f1*f9 + 2*(f3 + f7)^2 + 4*f1^2*foc1^2 + 2*f2^2*foc2^2 + 2*f4^2*foc2^2) + 4*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + f9*(2*f3 + 2*f7))*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f6 + f8)*(2*f3 + 2*f7) + f9*(2*f2 + 2*f4)) + (2*(f2 + f4)*(f3 + f7) + 2*f1*(f6 + f8))*(4*f7^2*foc1^2 + 4*f8^2*foc2^2 + 4*f9^2) + 8*foc2*((f1*f5 + f4*(f2 + f4))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + (f1*f6 + f4*(f3 + f7))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2)) + 8*foc1*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + 8*foc2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + 4*foc1^2*foc2^2*(2*(f1*f5 + f2*f4)*(f1*f6 + f3*f4) + 2*f1*f4*(f2*f6 + f3*f5)) + 24*f1^2*f2*f3*foc1^4 + 24*f4^2*f5*f6*foc2^4);
exp7_g= - ((f9*(2*f3 + 2*f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)^2 - 4*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2)^2 - 4*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2)^2 - 2*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + f9*(2*f3 + 2*f7))^2 - 2*foc1^4*(2*f1^2*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) + 4*f1^2*f3^2) - 2*foc2^4*(2*f4^2*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2) + 4*f4^2*f6^2) + (4*f1*f9 + 2*(f3 + f7)^2 + 4*f1^2*foc1^2 + 2*f2^2*foc2^2 + 2*f4^2*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (4*f7^2*foc1^2 + 4*f8^2*foc2^2 + 4*f9^2)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) - 8*foc1*(f1*f3 + f1*(f3 + f7))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) - 8*foc2*(f1*f6 + f4*(f3 + f7))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) - 4*foc1^2*foc2^2*((f1*f6 + f3*f4)^2 + 2*f1*f4*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6)));

exp8_a=(2*f5^3*(2*f2 + 2*f4));
exp8_b= (f5^2*(4*(f6 + f8)*(2*f2 + 2*f4) + 4*f5*(2*f3 + 2*f7)) - f5^2*(2*(f6 + f8)*(2*f2 + 2*f4) + 2*f5*(2*f3 + 2*f7)) + 4*f5^2*(f6 + f8)*(2*f2 + 2*f4));
exp8_c= (- f5^2*(2*(f6 + f8)*(2*f3 + 2*f7) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + f5^2*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f6 + f8)*(2*f3 + 2*f7) + 4*f9*(2*f2 + 2*f4)) + 8*f5^2*foc2^2*(f4*f5 + f5*(f2 + f4)) - 2*f5*(f6 + f8)*(2*(f6 + f8)*(2*f2 + 2*f4) + 2*f5*(2*f3 + 2*f7)) + 2*f5*(f6 + f8)*(4*(f6 + f8)*(2*f2 + 2*f4) + 4*f5*(2*f3 + 2*f7)) - 2*f5*(2*f2 + 2*f4)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 4*f5*(2*f2 + 2*f4)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 8*f2*f5*foc1^2*(f1*f5 + f2*(f2 + f4)));
exp8_d= (f5^2*(8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 4*f9*(2*f3 + 2*f7)) - (2*(f6 + f8)*(2*f2 + 2*f4) + 2*f5*(2*f3 + 2*f7))*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + (4*(f6 + f8)*(2*f2 + 2*f4) + 4*f5*(2*f3 + 2*f7))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) - f5^2*(2*f9*(2*f3 + 2*f7) + 4*f1*f3*foc1^2 + 4*f1*f7*foc1^2 + 4*f2*f8*foc2^2 + 4*f4*f6*foc2^2) + 8*foc1^2*(f1*f5 + f2*(f2 + f4))*(f3*f5 + f2*(f6 + f8)) + 8*foc2^2*(f4*f5 + f5*(f2 + f4))*(f5*f6 + f5*(f6 + f8)) + 8*f5^2*foc2^2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)) - 2*f5*(f6 + f8)*(2*(f6 + f8)*(2*f3 + 2*f7) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + 2*f5*(f6 + f8)*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f6 + f8)*(2*f3 + 2*f7) + 4*f9*(2*f2 + 2*f4)) - 2*f5*(2*f2 + 2*f4)*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 4*f5*(2*f2 + 2*f4)*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8)) + 8*f2*f5*foc1^2*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)));
exp8_e=  ((4*(f6 + f8)*(2*f2 + 2*f4) + 4*f5*(2*f3 + 2*f7))*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8)) - (2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2)*(2*(f6 + f8)*(2*f3 + 2*f7) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + (8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f6 + f8)*(2*f3 + 2*f7) + 4*f9*(2*f2 + 2*f4))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) - (2*(f6 + f8)*(2*f2 + 2*f4) + 2*f5*(2*f3 + 2*f7))*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 8*foc2*(f4*f5 + f5*(f2 + f4))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + 8*foc1*(f1*f5 + f2*(f2 + f4))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + 8*f1*f2^3*foc1^4 + 8*f4*f5^3*foc2^4 + f5^2*foc2*(8*foc2*(f4*f9 + f6*(f3 + f7)) + 8*f2*f5*foc2^3 + 8*f1*f4*foc1^2*foc2) + 4*f5*(2*f2 + 2*f4)*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 2*f5*(2*f2 + 2*f4)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 2*f5*(f6 + f8)*(2*f9*(2*f3 + 2*f7) + 4*f1*f3*foc1^2 + 4*f1*f7*foc1^2 + 4*f2*f8*foc2^2 + 4*f4*f6*foc2^2) + 8*foc1^2*(f3*f5 + f2*(f6 + f8))*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)) + 8*foc2^2*(f5*f6 + f5*(f6 + f8))*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)) + 2*f5*(f6 + f8)*(8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 4*f9*(2*f3 + 2*f7)) + f2*f5*foc1*(8*foc1*(f1*f9 + f3*(f3 + f7)) + 8*f1^2*foc1^3 + 8*f2^2*foc1*foc2^2) + 8*f2*f5*foc1^2*foc2^2*(f1*f5 + f2*f4));
exp8_f= (- (2*f9*(2*f3 + 2*f7) + 4*f1*f3*foc1^2 + 4*f1*f7*foc1^2 + 4*f2*f8*foc2^2 + 4*f4*f6*foc2^2)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + (8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 4*f9*(2*f3 + 2*f7))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) - (2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)*(2*(f6 + f8)*(2*f3 + 2*f7) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + (4*(f6 + f8)*(2*f2 + 2*f4) + 4*f5*(2*f3 + 2*f7))*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (2*(f6 + f8)*(2*f2 + 2*f4) + 2*f5*(2*f3 + 2*f7))*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8))*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f6 + f8)*(2*f3 + 2*f7) + 4*f9*(2*f2 + 2*f4)) + foc1*(f3*f5 + f2*(f6 + f8))*(8*foc1*(f1*f9 + f3*(f3 + f7)) + 8*f1^2*foc1^3 + 8*f2^2*foc1*foc2^2) + foc2*(f5*f6 + f5*(f6 + f8))*(8*foc2*(f4*f9 + f6*(f3 + f7)) + 8*f2*f5*foc2^3 + 8*f1*f4*foc1^2*foc2) + 8*foc1*(f1*f5 + f2*(f2 + f4))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + 8*foc2*(f4*f5 + f5*(f2 + f4))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + 8*foc2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + 8*foc1*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + 8*foc1^2*foc2^2*((f1*f5 + f2*f4)*(f2*f6 + f3*f5) + f2*f5*(f1*f6 + f3*f4)) + 24*f1*f2^2*f3*foc1^4 + 24*f4*f5^2*f6*foc2^4);
exp8_g= ((8*foc1*(f1*f9 + f3*(f3 + f7)) + 8*f1^2*foc1^3 + 8*f2^2*foc1*foc2^2)*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + (foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2)*(8*foc2*(f4*f9 + f6*(f3 + f7)) + 8*f2*f5*foc2^3 + 8*f1*f4*foc1^2*foc2) - (2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)*(2*f9*(2*f3 + 2*f7) + 4*f1*f3*foc1^2 + 4*f1*f7*foc1^2 + 4*f2*f8*foc2^2 + 4*f4*f6*foc2^2) + (2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8))*(8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 4*f9*(2*f3 + 2*f7)) + 4*foc1^4*(2*f1*f2*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) + 4*f1*f2*f3^2) + 4*foc2^4*(2*f4*f5*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2) + 4*f4*f5*f6^2) - (2*(f6 + f8)*(2*f3 + 2*f7) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (f7^2*foc1^2 + f8^2*foc2^2 + f9^2)*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f6 + f8)*(2*f3 + 2*f7) + 4*f9*(2*f2 + 2*f4)) + 8*foc1*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + 8*foc2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + 8*foc1^2*foc2^2*((f1*f5 + f2*f4)*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6) + (f1*f6 + f3*f4)*(f2*f6 + f3*f5)));
exp8_h= - ((2*f9*(2*f3 + 2*f7) + 4*f1*f3*foc1^2 + 4*f1*f7*foc1^2 + 4*f2*f8*foc2^2 + 4*f4*f6*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 4*f9*(2*f3 + 2*f7))*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1)*(8*foc1*(f1*f9 + f3*(f3 + f7)) + 8*f1^2*foc1^3 + 8*f2^2*foc1*foc2^2) - (f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2)*(8*foc2*(f4*f9 + f6*(f3 + f7)) + 8*f2*f5*foc2^3 + 8*f1*f4*foc1^2*foc2) - 8*foc1^2*foc2^2*(f1*f6 + f3*f4)*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6) - 8*f1*f3*foc1^4*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) - 8*f4*f6*foc2^4*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2)) ;

exp9_a=f5^4;
exp9_b=(4*f5^3*(f6 + f8));
exp9_c=  (4*f5^2*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 4*f5^2*(f6 + f8)^2 + 4*f5^4*foc2^2 - 2*f5^2*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 4*f2^2*f5^2*foc1^2);
exp9_d= (- 2*f5^2*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 4*f5^2*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8)) + 8*f5^2*foc2^2*(f5*f6 + f5*(f6 + f8)) - 4*f5*(f6 + f8)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 8*f5*(f6 + f8)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 8*f2*f5*foc1^2*(f3*f5 + f2*(f6 + f8)));
exp9_e= (4*foc1^2*(f3*f5 + f2*(f6 + f8))^2 + 4*foc2^2*(f5*f6 + f5*(f6 + f8))^2 + 4*f5^2*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 2*f5^2*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2)^2 + 2*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2)^2 + 2*f2^4*foc1^4 + 2*f5^4*foc2^4 - 4*f5*(f6 + f8)*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 8*f5*(f6 + f8)*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8)) + 8*f5^2*foc2*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + 4*f2^2*f5^2*foc1^2*foc2^2 + 8*f2*f5*foc1*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1));
exp9_f=  (- 2*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 4*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 8*foc2*(f5*f6 + f5*(f6 + f8))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + 8*foc1*(f3*f5 + f2*(f6 + f8))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + 8*f2^3*f3*foc1^4 + 8*f5^3*f6*foc2^4 + 8*f5^2*foc2*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + 8*f5*(f6 + f8)*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 4*f5*(f6 + f8)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 8*f2*f5*foc1*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + 8*f2*f5*foc1^2*foc2^2*(f2*f6 + f3*f5));
exp9_g= (4*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2)^2 - 2*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 4*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 4*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1)^2 - (2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)^2 + 2*foc1^4*(2*f2^2*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) + 4*f2^2*f3^2) + 2*foc2^4*(2*f5^2*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2) + 4*f5^2*f6^2) + 2*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8))^2 + 8*foc1*(f3*f5 + f2*(f6 + f8))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + 8*foc2*(f5*f6 + f5*(f6 + f8))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + 4*foc1^2*foc2^2*((f2*f6 + f3*f5)^2 + 2*f2*f5*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6)));
exp9_h=  (8*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2)*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) - 2*(2*f9*(f6 + f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 8*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1)*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + 4*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + 2*f9*(f6 + f8))*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 8*foc1^2*foc2^2*(f2*f6 + f3*f5)*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6) + 8*f2*f3*foc1^4*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) + 8*f5*f6*foc2^4*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2));
exp9_i= (2*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2)^2 - (f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2)^2 + 2*foc1^4*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2)^2 + 2*foc2^4*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2)^2 + 4*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1)^2 + 4*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2)^2 + 4*foc1^2*foc2^2*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6)^2);
for i=1:numsols
    
    yc=xarray(i);
    
    yc2=yc^2;
    yc3=yc*yc2;
    yc4=yc*yc3 ;
    yc5=yc*yc4 ;
    yc6=yc*yc5 ;
    yc7=yc*yc6 ;
    yc8=yc*yc7 ;
    
    a1= exp1       ;
    b1=(exp2_a*yc )+ exp2_b ;
    c1=(exp3_a*yc2 )+ (exp3_b*yc )+ exp3_c;
    d1=(exp4_a*yc3 )+(exp4_b*yc2 )+ (exp4_c*yc  )+ exp4_d ;
    e1=(exp5_a*yc4 )+(exp5_b*yc3 )+ (exp5_c*yc2 )+   (exp5_d*yc ) + exp5_e ;
    f1=(exp6_a*yc5 )+(exp6_b*yc4 )+ (exp6_c*yc3 )+   (exp6_d*yc2 ) +(exp6_e*yc )+   exp6_f ;
    g1=(exp7_a*yc6 )+(exp7_b*yc5 )+ (exp7_c*yc4 )+   (exp7_d*yc3 ) +(exp7_e*yc2 )+   (exp7_f*yc )+    exp7_g  ;
    h1=(exp8_a*yc7 )+(exp8_b*yc6 )+ (exp8_c*yc5 )+   (exp8_d*yc4 ) +(exp8_e*yc3 )+   (exp8_f*yc2 )+   (exp8_g*yc )+ exp8_h  ;
    i1=(exp9_a*yc8 )+(exp9_b*yc7 )+ (exp9_c*yc6 )+   (exp9_d*yc5 ) +(exp9_e*yc4 )+   (exp9_f*yc3 )+   (exp9_g*yc2 )+  (exp9_h*yc )+ exp9_i  ;
    
    p=[a1; b1 ; c1; d1; e1; f1; g1; h1 ; i1    ];
    
  
  f = @(x)((p(1)*x.^8)+(p(2)*x.^7)+(p(3)*x.^6)+(p(4)*x.^5)+(p(5)*x.^4)+(p(6)*x.^3)+(p(7)*x.^2)+(p(8)*x.^1)+(p(9)*x.^0));
  answ=processRoots_focalonly(fminbnd(f, xcstart, xcend));
    if(answ>minans && answ<maxans)
        solzarray=[solzarray ;answ];
        xfarray= [xfarray ; xarray(i)];
        
    end
    
end
xarray=xarray';


end



function [solzarray,xarray,xfarray]=findOCyc(F,w,h,numsols,Kguess)

F=F/norm(F,'fro');

ycstart=10;
ycend=h;
if(numsols==1)
    xarray=[Kguess(2,2)];
else
    step=((ycend-ycstart)/(numsols-1));
    xarray=ycstart:step:ycend;
    
end

solzarray=[];
xfarray=[];




minans=1;
maxans=ycend;

f1=F(1,1);
f2=F(1,2);
f3=F(1,3);

f4=F(2,1);
f5=F(2,2);
f6=F(2,3);

f7=F(3,1);
f8=F(3,2);
f9=F(3,3);


foc1=Kguess(1,1);
foc2=Kguess(2,2);



exp1= f5^4;

exp2_a= (4*f5^3*(f2 + f4));
exp2_b= 4*f5^3*(f6 + f8);

exp3_a= (2*f5^2*(2*f1*f5 + (f2 + f4)^2) + 4*f5^2*(f2 + f4)^2);
exp3_b= (2*f5^2*(2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7)) + 8*f5^2*(f2 + f4)*(f6 + f8));
exp3_c= 4*f5^2*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 4*f5^2*(f6 + f8)^2 + 4*f5^4*foc2^2 - 2*f5^2*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 4*f2^2*f5^2*foc1^2;


exp4_a= (4*f5*(2*f1*f5 + (f2 + f4)^2)*(f2 + f4) + 2*f1*f5^2*(2*f2 + 2*f4)) ;
exp4_b= (4*f5*((2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(f2 + f4) + (2*f1*f5 + (f2 + f4)^2)*(f6 + f8)) + 2*f5^2*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8))) ;
exp4_c= (- 4*f5*((f2 + f4)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(f6 + f8)) + 8*f5*((f2 + f4)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(f6 + f8)) - 2*f5^2*((f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2) + 4*f5^2*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4)) + 8*f5^2*foc2^2*(f4*f5 + f5*(f2 + f4)) + 8*f2*f5*foc1^2*(f1*f5 + f2*(f2 + f4)));
exp4_d=  4*f5^2*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + f9*(2*f6 + 2*f8)) - 2*f5^2*(f9*(2*f6 + 2*f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 8*f5^2*foc2^2*(f5*f6 + f5*(f6 + f8)) - 4*f5*(f6 + f8)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 8*f5*(f6 + f8)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 8*f2*f5*foc1^2*(f3*f5 + f2*(f6 + f8));

exp5_a=((2*f1*f5 + (f2 + f4)^2)^2 + 2*f1^2*f5^2 + 4*f1*f5*(f2 + f4)*(2*f2 + 2*f4));
exp5_b=  (4*f5*((f2 + f4)*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8)) + f1*(f6 + f8)*(2*f2 + 2*f4)) + 2*(2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(2*f1*f5 + (f2 + f4)^2) + 4*f1*f5^2*(f3 + f7));
exp5_c=  (4*foc1^2*(f1*f5 + f2*(f2 + f4))^2 + 4*foc2^2*(f4*f5 + f5*(f2 + f4))^2 - 4*f5*((f6 + f8)*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8)) + (f2 + f4)*((f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2)) + 4*f5^2*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 8*f5*((f6 + f8)*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8)) + (f2 + f4)*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4))) - 2*(2*f1*f5 + (f2 + f4)^2)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 4*(2*f1*f5 + (f2 + f4)^2)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))^2 - 2*f5^2*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 8*f5^2*foc2^2*(f1*f5 + f4*(f2 + f4)) + 8*f2*f5*foc1^2*(f1*f2 + f1*(f2 + f4)));
exp5_d=  (- 2*f5^2*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) - 2*(2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 4*(2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 4*f5^2*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7)) - 4*f5*((f2 + f4)*(f9*(2*f6 + 2*f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + (f6 + f8)*((f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2)) + 8*f5*((f2 + f4)*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + f9*(2*f6 + 2*f8)) + (f6 + f8)*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4))) + 8*foc1^2*(f1*f5 + f2*(f2 + f4))*(f3*f5 + f2*(f6 + f8)) + 8*foc2^2*(f4*f5 + f5*(f2 + f4))*(f5*f6 + f5*(f6 + f8)) + 8*f5^2*foc2^2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)) + 8*f2*f5*foc1^2*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)));
exp5_e= (4*foc1^2*(f3*f5 + f2*(f6 + f8))^2 + 4*foc2^2*(f5*f6 + f5*(f6 + f8))^2 + 4*f5^2*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 2*f5^2*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2)^2 + 2*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2)^2 + 2*f2^4*foc1^4 + 2*f5^4*foc2^4 + 8*f5^2*foc2*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) - 4*f5*(f6 + f8)*(f9*(2*f6 + 2*f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 8*f5*(f6 + f8)*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + f9*(2*f6 + 2*f8)) + 4*f2^2*f5^2*foc1^2*foc2^2 + 8*f2*f5*foc1*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1));

exp6_a= (4*f1^2*f5*(f2 + f4) + 2*f1*(2*f1*f5 + (f2 + f4)^2)*(2*f2 + 2*f4));
exp6_b=(4*f5*((f6 + f8)*f1^2 + 2*(f2 + f4)*(f3 + f7)*f1) - (2*f1*f5 + (f2 + f4)^2)*(2*(f3 + f7)*(2*f2 + 2*f4) + 2*f1*(2*f6 + 2*f8)) + (2*f1*f5 + (f2 + f4)^2)*(4*(f3 + f7)*(2*f2 + 2*f4) + 4*f1*(2*f6 + 2*f8)) + 2*f1*(2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(2*f2 + 2*f4));
exp6_c= (- 4*f5*((f2 + f4)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 2*f1*(f3 + f7)*(f6 + f8)) + 8*f5*((f2 + f4)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 2*f1*(f3 + f7)*(f6 + f8)) - (2*f1*f5 + (f2 + f4)^2)*(2*(f3 + f7)*(2*f6 + 2*f8) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + (2*f1*f5 + (f2 + f4)^2)*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f3 + f7)*(2*f6 + 2*f8) + 4*f9*(2*f2 + 2*f4)) - (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(2*(f3 + f7)*(2*f2 + 2*f4) + 2*f1*(2*f6 + 2*f8)) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(4*(f3 + f7)*(2*f2 + 2*f4) + 4*f1*(2*f6 + 2*f8)) + 8*foc1^2*(f1*f2 + f1*(f2 + f4))*(f1*f5 + f2*(f2 + f4)) + 8*foc2^2*(f1*f5 + f4*(f2 + f4))*(f4*f5 + f5*(f2 + f4)) - 2*f1*(2*f2 + 2*f4)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + 4*f1*(2*f2 + 2*f4)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 8*f1^2*f2*f5*foc1^2 + 8*f1*f4*f5^2*foc2^2);
exp6_d=  (- 4*f5*((f6 + f8)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + (f2 + f4)*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)) - (2*f1*f5 + (f2 + f4)^2)*(2*f9*(2*f6 + 2*f8) + 4*f2*f3*foc1^2 + 4*f4*f7*foc1^2 + 4*f5*f6*foc2^2 + 4*f5*f8*foc2^2) + 8*foc1*(foc1*(f1*f5 + f2*(f2 + f4))*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)) + foc1*(f1*f2 + f1*(f2 + f4))*(f3*f5 + f2*(f6 + f8))) + 8*foc2*(foc2*(f4*f5 + f5*(f2 + f4))*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)) + foc2*(f1*f5 + f4*(f2 + f4))*(f5*f6 + f5*(f6 + f8))) + 8*f5*((f6 + f8)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + (f2 + f4)*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7))) + (2*f1*f5 + (f2 + f4)^2)*(8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 4*f9*(2*f6 + 2*f8)) - (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(2*(f3 + f7)*(2*f6 + 2*f8) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) - (2*(f3 + f7)*(2*f2 + 2*f4) + 2*f1*(2*f6 + 2*f8))*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f3 + f7)*(2*f6 + 2*f8) + 4*f9*(2*f2 + 2*f4)) + (4*(f3 + f7)*(2*f2 + 2*f4) + 4*f1*(2*f6 + 2*f8))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 8*f5^2*foc2^2*(f1*f6 + f4*(f3 + f7)) + 8*f2*f5*foc1^2*(f1*f3 + f1*(f3 + f7)));
exp6_e=  (- (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(2*f9*(2*f6 + 2*f8) + 4*f2*f3*foc1^2 + 4*f4*f7*foc1^2 + 4*f5*f6*foc2^2 + 4*f5*f8*foc2^2) + 8*foc2*((f4*f5 + f5*(f2 + f4))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + foc2*(f5*f6 + f5*(f6 + f8))*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 4*f9*(2*f6 + 2*f8)) - 4*f5*((f2 + f4)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (f6 + f8)*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)) + 8*foc1*((f1*f5 + f2*(f2 + f4))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + foc1*(f3*f5 + f2*(f6 + f8))*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))) + 8*f5*((f2 + f4)*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (f6 + f8)*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7))) - (2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2)*(2*(f3 + f7)*(2*f6 + 2*f8) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + (8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f3 + f7)*(2*f6 + 2*f8) + 4*f9*(2*f2 + 2*f4))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + 8*f1*f2^3*foc1^4 + 8*f4*f5^3*foc2^4 + 8*f5^2*foc2*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + 8*f2*f5*foc1*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + 8*f2*f5*foc1^2*foc2^2*(f1*f5 + f2*f4));
exp6_f= - ((2*f9*(2*f6 + 2*f8) + 4*f2*f3*foc1^2 + 4*f4*f7*foc1^2 + 4*f5*f6*foc2^2 + 4*f5*f8*foc2^2)*(2*f5*f9 + (f6 + f8)^2 + f2^2*foc1^2 + f4^2*foc1^2 + 2*f5^2*foc2^2) - (8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 4*f9*(2*f6 + 2*f8))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) - 8*foc2*(f5*f6 + f5*(f6 + f8))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) - 8*foc1*(f3*f5 + f2*(f6 + f8))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) - 8*f2^3*f3*foc1^4 - 8*f5^3*f6*foc2^4 - 8*f5^2*foc2*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) - 8*f5*(f6 + f8)*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 4*f5*(f6 + f8)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 8*f2*f5*foc1*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) - 8*f2*f5*foc1^2*foc2^2*(f2*f6 + f3*f5));

exp7_a= (4*f1^2*(2*f1*f5 + (f2 + f4)^2) + f1^2*(2*f2 + 2*f4)^2 - f1^2*(4*f1*f5 + 2*(f2 + f4)^2));
exp7_b=(4*f1^2*(2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7)) - f1^2*(4*(f2 + f4)*(f6 + f8) + 4*f5*(f3 + f7)) + 8*f1*(2*f1*f5 + (f2 + f4)^2)*(f3 + f7) + 2*f1*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8))*(2*f2 + 2*f4) - 2*f1*(f3 + f7)*(4*f1*f5 + 2*(f2 + f4)^2));
exp7_c=  (4*foc1^2*(f1*f2 + f1*(f2 + f4))^2 + 4*foc2^2*(f1*f5 + f4*(f2 + f4))^2 + 4*f1^2*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) - f1^2*(4*f5*f9 + 2*(f6 + f8)^2 + 2*f2^2*foc1^2 + 2*f4^2*foc1^2 + 4*f5^2*foc2^2) - (4*f1*f5 + 2*(f2 + f4)^2)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + (2*f1*f5 + (f2 + f4)^2)*(8*f1*f9 + 4*(f3 + f7)^2 + 4*f1^2*foc1^2 + 4*f2^2*foc2^2) + ((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8))^2 + 4*f1*(2*f2 + 2*f4)*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4)) + 8*f1^2*foc1^2*(f1*f5 + f2*(f2 + f4)) + 8*f1*(2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(f3 + f7) - 2*f1*(4*(f2 + f4)*(f6 + f8) + 4*f5*(f3 + f7))*(f3 + f7) - 2*f1*(2*f2 + 2*f4)*((f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2) + 8*f1*f4*foc2^2*(f4*f5 + f5*(f2 + f4)));
exp7_d= (- 2*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8))*((f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2) - (4*(f2 + f4)*(f6 + f8) + 4*f5*(f3 + f7))*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 4*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8))*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4)) + (2*f1*f5 + (f2 + f4)^2)*(8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 8*f9*(f3 + f7)) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(8*f1*f9 + 4*(f3 + f7)^2 + 4*f1^2*foc1^2 + 4*f2^2*foc2^2) + 8*foc1*(f1^2*foc1*(f3*f5 + f2*(f6 + f8)) + foc1*(f1*f5 + f2*(f2 + f4))*(f1*f3 + f1*(f3 + f7))) - (4*f1*f5 + 2*(f2 + f4)^2)*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 8*foc2*(foc2*(f4*f5 + f5*(f2 + f4))*(f1*f6 + f4*(f3 + f7)) + f1*f4*foc2*(f5*f6 + f5*(f6 + f8))) - 2*f1*(f3 + f7)*(4*f5*f9 + 2*(f6 + f8)^2 + 2*f2^2*foc1^2 + 2*f4^2*foc1^2 + 4*f5^2*foc2^2) - 2*f1*(2*f2 + 2*f4)*(f9*(2*f6 + 2*f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 4*f1*(2*f2 + 2*f4)*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + f9*(2*f6 + 2*f8)) + 8*foc1^2*(f1*f2 + f1*(f2 + f4))*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)) + 8*foc2^2*(f1*f5 + f4*(f2 + f4))*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)) + 8*f1*(f3 + f7)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2));
exp7_e= (- (4*f5*f9 + 2*(f6 + f8)^2 + 2*f2^2*foc1^2 + 2*f4^2*foc1^2 + 4*f5^2*foc2^2)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) - (4*(f2 + f4)*(f6 + f8) + 4*f5*(f3 + f7))*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + (8*f1*f9 + 4*(f3 + f7)^2 + 4*f1^2*foc1^2 + 4*f2^2*foc2^2)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 8*f9*(f3 + f7)) + 8*foc1*((f1*f5 + f2*(f2 + f4))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + foc1*(f1*f3 + f1*(f3 + f7))*(f3*f5 + f2*(f6 + f8))) - ((f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2)^2 + 8*foc2*((f4*f5 + f5*(f2 + f4))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + foc2*(f1*f6 + f4*(f3 + f7))*(f5*f6 + f5*(f6 + f8))) + 4*foc1^2*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))^2 + 4*foc2^2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))^2 - 2*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8))*(f9*(2*f6 + 2*f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2) + 2*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4))^2 + (2*f1*f5 + (f2 + f4)^2)*(4*f7^2*foc1^2 + 4*f8^2*foc2^2 + 4*f9^2) + 4*((f3 + f7)*(2*f2 + 2*f4) + f1*(2*f6 + 2*f8))*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + f9*(2*f6 + 2*f8)) - (4*f1*f5 + 2*(f2 + f4)^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 8*foc2*(f1*f5 + f4*(f2 + f4))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + 8*foc1*(f1*f2 + f1*(f2 + f4))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + 4*foc1^2*foc2^2*((f1*f5 + f2*f4)^2 + 2*f1*f2*f4*f5) + 12*f1^2*f2^2*foc1^4 + 12*f4^2*f5^2*foc2^4);
exp7_f=  (8*foc2*((f4*f5 + f5*(f2 + f4))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + (f5*f6 + f5*(f6 + f8))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2)) - (4*(f2 + f4)*(f6 + f8) + 4*f5*(f3 + f7))*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (8*f1*f7*foc1^2 + 8*f2*f8*foc2^2 + 8*f9*(f3 + f7))*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) - 2*(f9*(2*f6 + 2*f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)*((f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4) + 2*f1*f2*foc1^2 + 2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + 2*f4*f5*foc2^2) - (2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)*(4*f5*f9 + 2*(f6 + f8)^2 + 2*f2^2*foc1^2 + 2*f4^2*foc1^2 + 4*f5^2*foc2^2) + 4*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + f9*(2*f6 + 2*f8))*(2*f1*f4*foc1^2 + 2*f2*f5*foc2^2 + (f3 + f7)*(2*f6 + 2*f8) + f9*(2*f2 + 2*f4)) + (2*(f2 + f4)*(f6 + f8) + 2*f5*(f3 + f7))*(4*f7^2*foc1^2 + 4*f8^2*foc2^2 + 4*f9^2) + 8*foc1*((f1*f5 + f2*(f2 + f4))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + (f3*f5 + f2*(f6 + f8))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2)) + 8*foc2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2) + 8*foc1*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1) + 4*foc1^2*foc2^2*(2*(f1*f5 + f2*f4)*(f2*f6 + f3*f5) + 2*f2*f5*(f1*f6 + f3*f4)) + 24*f1*f2^2*f3*foc1^4 + 24*f4*f5^2*f6*foc2^4);
exp7_g=  - ((f9*(2*f6 + 2*f8) + 2*f2*f3*foc1^2 + 2*f4*f7*foc1^2 + 2*f5*f6*foc2^2 + 2*f5*f8*foc2^2)^2 - 4*(foc2*(f5*f9 + f6*(f6 + f8)) + f5^2*foc2^3 + f4^2*foc1^2*foc2)^2 - 4*(f1*f4*foc1^3 + f2*f5*foc1*foc2^2 + (f2*f9 + f3*(f6 + f8))*foc1)^2 - 2*(2*f4*f7*foc1^2 + 2*f5*f8*foc2^2 + f9*(2*f6 + 2*f8))^2 - 2*foc1^4*(2*f2^2*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) + 4*f2^2*f3^2) - 2*foc2^4*(2*f5^2*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2) + 4*f5^2*f6^2) + (4*f5*f9 + 2*(f6 + f8)^2 + 2*f2^2*foc1^2 + 2*f4^2*foc1^2 + 4*f5^2*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (4*f7^2*foc1^2 + 4*f8^2*foc2^2 + 4*f9^2)*(2*f5*f9 + (f6 + f8)^2 + f4^2*foc1^2 + f5^2*foc2^2) - 8*foc1*(f3*f5 + f2*(f6 + f8))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) - 8*foc2*(f5*f6 + f5*(f6 + f8))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) - 4*foc1^2*foc2^2*((f2*f6 + f3*f5)^2 + 2*f2*f5*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6)));

exp8_a=(2*f1^3*(2*f2 + 2*f4)) ;
exp8_b=(f1^2*(4*(f3 + f7)*(2*f2 + 2*f4) + 4*f1*(2*f6 + 2*f8)) - f1^2*(2*(f3 + f7)*(2*f2 + 2*f4) + 2*f1*(2*f6 + 2*f8)) + 4*f1^2*(f3 + f7)*(2*f2 + 2*f4)) ;
exp8_c=  (- f1^2*(2*(f3 + f7)*(2*f6 + 2*f8) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + f1^2*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f3 + f7)*(2*f6 + 2*f8) + 4*f9*(2*f2 + 2*f4)) + 8*f1^2*foc1^2*(f1*f2 + f1*(f2 + f4)) - 2*f1*(f3 + f7)*(2*(f3 + f7)*(2*f2 + 2*f4) + 2*f1*(2*f6 + 2*f8)) + 2*f1*(f3 + f7)*(4*(f3 + f7)*(2*f2 + 2*f4) + 4*f1*(2*f6 + 2*f8)) - 2*f1*(2*f2 + 2*f4)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 4*f1*(2*f2 + 2*f4)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 8*f1*f4*foc2^2*(f1*f5 + f4*(f2 + f4)));
exp8_d= (f1^2*(8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 4*f9*(2*f6 + 2*f8)) - (2*(f3 + f7)*(2*f2 + 2*f4) + 2*f1*(2*f6 + 2*f8))*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + (4*(f3 + f7)*(2*f2 + 2*f4) + 4*f1*(2*f6 + 2*f8))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) - f1^2*(2*f9*(2*f6 + 2*f8) + 4*f2*f3*foc1^2 + 4*f4*f7*foc1^2 + 4*f5*f6*foc2^2 + 4*f5*f8*foc2^2) + 8*foc1^2*(f1*f2 + f1*(f2 + f4))*(f1*f3 + f1*(f3 + f7)) + 8*foc2^2*(f1*f5 + f4*(f2 + f4))*(f1*f6 + f4*(f3 + f7)) + 8*f1^2*foc1^2*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)) - 2*f1*(f3 + f7)*(2*(f3 + f7)*(2*f6 + 2*f8) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + 2*f1*(f3 + f7)*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f3 + f7)*(2*f6 + 2*f8) + 4*f9*(2*f2 + 2*f4)) - 2*f1*(2*f2 + 2*f4)*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 4*f1*(2*f2 + 2*f4)*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7)) + 8*f1*f4*foc2^2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)));
exp8_e= ((4*(f3 + f7)*(2*f2 + 2*f4) + 4*f1*(2*f6 + 2*f8))*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7)) - (2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2)*(2*(f3 + f7)*(2*f6 + 2*f8) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + (8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f3 + f7)*(2*f6 + 2*f8) + 4*f9*(2*f2 + 2*f4))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) - (2*(f3 + f7)*(2*f2 + 2*f4) + 2*f1*(2*f6 + 2*f8))*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 8*foc1*(f1*f2 + f1*(f2 + f4))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + 8*foc2*(f1*f5 + f4*(f2 + f4))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + 8*f1^3*f2*foc1^4 + 8*f4^3*f5*foc2^4 + f1^2*foc1*(8*foc1*(f2*f9 + f3*(f6 + f8)) + 8*f1*f4*foc1^3 + 8*f2*f5*foc1*foc2^2) + 4*f1*(2*f2 + 2*f4)*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 2*f1*(2*f2 + 2*f4)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 2*f1*(f3 + f7)*(2*f9*(2*f6 + 2*f8) + 4*f2*f3*foc1^2 + 4*f4*f7*foc1^2 + 4*f5*f6*foc2^2 + 4*f5*f8*foc2^2) + 8*foc1^2*(f1*f3 + f1*(f3 + f7))*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8)) + 8*foc2^2*(f1*f6 + f4*(f3 + f7))*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8)) + 2*f1*(f3 + f7)*(8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 4*f9*(2*f6 + 2*f8)) + f1*f4*foc2*(8*foc2*(f5*f9 + f6*(f6 + f8)) + 8*f5^2*foc2^3 + 8*f4^2*foc1^2*foc2) + 8*f1*f4*foc1^2*foc2^2*(f1*f5 + f2*f4));
exp8_f=  (- (2*f9*(2*f6 + 2*f8) + 4*f2*f3*foc1^2 + 4*f4*f7*foc1^2 + 4*f5*f6*foc2^2 + 4*f5*f8*foc2^2)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + (8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 4*f9*(2*f6 + 2*f8))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) - (2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)*(2*(f3 + f7)*(2*f6 + 2*f8) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2) + (4*(f3 + f7)*(2*f2 + 2*f4) + 4*f1*(2*f6 + 2*f8))*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (2*(f3 + f7)*(2*f2 + 2*f4) + 2*f1*(2*f6 + 2*f8))*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7))*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f3 + f7)*(2*f6 + 2*f8) + 4*f9*(2*f2 + 2*f4)) + foc2*(f1*f6 + f4*(f3 + f7))*(8*foc2*(f5*f9 + f6*(f6 + f8)) + 8*f5^2*foc2^3 + 8*f4^2*foc1^2*foc2) + foc1*(f1*f3 + f1*(f3 + f7))*(8*foc1*(f2*f9 + f3*(f6 + f8)) + 8*f1*f4*foc1^3 + 8*f2*f5*foc1*foc2^2) + 8*foc1*(f1*f2 + f1*(f2 + f4))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + 8*foc2*(f1*f5 + f4*(f2 + f4))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + 8*foc1*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + 8*foc2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + 8*foc1^2*foc2^2*((f1*f5 + f2*f4)*(f1*f6 + f3*f4) + f1*f4*(f2*f6 + f3*f5)) + 24*f1^2*f2*f3*foc1^4 + 24*f4^2*f5*f6*foc2^4);
exp8_g= ((foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2)*(8*foc1*(f2*f9 + f3*(f6 + f8)) + 8*f1*f4*foc1^3 + 8*f2*f5*foc1*foc2^2) + (8*foc2*(f5*f9 + f6*(f6 + f8)) + 8*f5^2*foc2^3 + 8*f4^2*foc1^2*foc2)*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) - (2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)*(2*f9*(2*f6 + 2*f8) + 4*f2*f3*foc1^2 + 4*f4*f7*foc1^2 + 4*f5*f6*foc2^2 + 4*f5*f8*foc2^2) + (2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7))*(8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 4*f9*(2*f6 + 2*f8)) + 4*foc1^4*(2*f1*f2*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) + 4*f1*f2*f3^2) + 4*foc2^4*(2*f4*f5*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2) + 4*f4*f5*f6^2) - (2*(f3 + f7)*(2*f6 + 2*f8) + 2*f9*(2*f2 + 2*f4) + 4*f1*f2*foc1^2 + 4*f1*f4*foc1^2 + 4*f2*f5*foc2^2 + 4*f4*f5*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + (f7^2*foc1^2 + f8^2*foc2^2 + f9^2)*(8*f1*f4*foc1^2 + 8*f2*f5*foc2^2 + 4*(f3 + f7)*(2*f6 + 2*f8) + 4*f9*(2*f2 + 2*f4)) + 8*foc1*(f3*(f2 + f4) + f2*(f3 + f7) + f1*(f6 + f8))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + 8*foc2*(f6*(f2 + f4) + f5*(f3 + f7) + f4*(f6 + f8))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + 8*foc1^2*foc2^2*((f1*f5 + f2*f4)*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6) + (f1*f6 + f3*f4)*(f2*f6 + f3*f5))) ;
exp8_h= - ((2*f9*(2*f6 + 2*f8) + 4*f2*f3*foc1^2 + 4*f4*f7*foc1^2 + 4*f5*f6*foc2^2 + 4*f5*f8*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (8*f4*f7*foc1^2 + 8*f5*f8*foc2^2 + 4*f9*(2*f6 + 2*f8))*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2)*(8*foc2*(f5*f9 + f6*(f6 + f8)) + 8*f5^2*foc2^3 + 8*f4^2*foc1^2*foc2) - (f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1)*(8*foc1*(f2*f9 + f3*(f6 + f8)) + 8*f1*f4*foc1^3 + 8*f2*f5*foc1*foc2^2) - 8*foc1^2*foc2^2*(f2*f6 + f3*f5)*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6) - 8*f2*f3*foc1^4*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) - 8*f5*f6*foc2^4*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2));

exp9_a= f1^4;
exp9_b=(4*f1^3*(f3 + f7));
exp9_c= (4*f1^2*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 4*f1^2*(f3 + f7)^2 + 4*f1^4*foc1^2 - 2*f1^2*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 4*f1^2*f4^2*foc2^2);
exp9_d= (- 2*f1^2*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 4*f1^2*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7)) + 8*f1^2*foc1^2*(f1*f3 + f1*(f3 + f7)) - 4*f1*(f3 + f7)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 8*f1*(f3 + f7)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 8*f1*f4*foc2^2*(f1*f6 + f4*(f3 + f7)));
exp9_e=  (4*foc1^2*(f1*f3 + f1*(f3 + f7))^2 + 4*foc2^2*(f1*f6 + f4*(f3 + f7))^2 + 4*f1^2*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 2*f1^2*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - (2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2)^2 + 2*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2)^2 + 2*f1^4*foc1^4 + 2*f4^4*foc2^4 - 4*f1*(f3 + f7)*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2) + 8*f1*(f3 + f7)*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7)) + 8*f1^2*foc1*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + 4*f1^2*f4^2*foc1^2*foc2^2 + 8*f1*f4*foc2*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2));
exp9_f=  (- 2*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2) + 4*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7))*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 8*foc1*(f1*f3 + f1*(f3 + f7))*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) + 8*foc2*(f1*f6 + f4*(f3 + f7))*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + 8*f1^3*f3*foc1^4 + 8*f4^3*f6*foc2^4 + 8*f1^2*foc1*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + 8*f1*(f3 + f7)*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) - 4*f1*(f3 + f7)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 8*f1*f4*foc2*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + 8*f1*f4*foc1^2*foc2^2*(f1*f6 + f3*f4));
exp9_g=  (4*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2)^2 - 2*(2*f1*f9 + (f3 + f7)^2 + 2*f1^2*foc1^2 + f2^2*foc2^2 + f4^2*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 4*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2)*(2*f1*f9 + (f3 + f7)^2 + f1^2*foc1^2 + f2^2*foc2^2) + 4*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2)^2 - (2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)^2 + 2*foc1^4*(2*f1^2*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) + 4*f1^2*f3^2) + 2*foc2^4*(2*f4^2*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2) + 4*f4^2*f6^2) + 2*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7))^2 + 8*foc1*(f1*f3 + f1*(f3 + f7))*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1) + 8*foc2*(f1*f6 + f4*(f3 + f7))*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2) + 4*foc1^2*foc2^2*((f1*f6 + f3*f4)^2 + 2*f1*f4*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6)));
exp9_h= (8*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1)*(foc1*(f1*f9 + f3*(f3 + f7)) + f1^2*foc1^3 + f2^2*foc1*foc2^2) - 2*(2*f9*(f3 + f7) + 2*f1*f3*foc1^2 + 2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f4*f6*foc2^2)*(f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 8*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2)*(f1*f4*foc1^2*foc2 + f2*f5*foc2^3 + (f4*f9 + f6*(f3 + f7))*foc2) + 4*(2*f1*f7*foc1^2 + 2*f2*f8*foc2^2 + 2*f9*(f3 + f7))*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2) + 8*foc1^2*foc2^2*(f1*f6 + f3*f4)*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6) + 8*f1*f3*foc1^4*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2) + 8*f4*f6*foc2^4*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2)) ;
exp9_i= (2*(f7^2*foc1^2 + f8^2*foc2^2 + f9^2)^2 - (f1^2*foc1^4 + f2^2*foc1^2*foc2^2 + f3^2*foc1^2 + f4^2*foc1^2*foc2^2 + f5^2*foc2^4 + f6^2*foc2^2 + f7^2*foc1^2 + f8^2*foc2^2 + f9^2)^2 + 2*foc1^4*(f1^2*foc1^2 + f2^2*foc2^2 + f3^2)^2 + 2*foc2^4*(f4^2*foc1^2 + f5^2*foc2^2 + f6^2)^2 + 4*(f1*f7*foc1^3 + f2*f8*foc1*foc2^2 + f3*f9*foc1)^2 + 4*(f4*f7*foc1^2*foc2 + f5*f8*foc2^3 + f6*f9*foc2)^2 + 4*foc1^2*foc2^2*(f1*f4*foc1^2 + f2*f5*foc2^2 + f3*f6)^2);

for i=1:numsols
    
    xc=xarray(i);
    xc2=xc^2;
    xc3=xc*xc2;
    xc4=xc*xc3 ;
    xc5=xc*xc4 ;
    xc6=xc*xc5 ;
    xc7=xc*xc6 ;
    xc8=xc*xc7 ;
    
    a1= exp1       ;
    b1=(exp2_a*xc )+ exp2_b ;
    c1=(exp3_a*xc2 )+ (exp3_b*xc )+ exp3_c;
    d1=(exp4_a*xc3 )+(exp4_b*xc2 )+ (exp4_c*xc  )+ exp4_d ;
    e1=(exp5_a*xc4 )+(exp5_b*xc3 )+ (exp5_c*xc2 )+   (exp5_d*xc ) + exp5_e ;
    f1=(exp6_a*xc5 )+(exp6_b*xc4 )+ (exp6_c*xc3 )+   (exp6_d*xc2 ) +(exp6_e*xc )+   exp6_f ;
    g1=(exp7_a*xc6 )+(exp7_b*xc5 )+ (exp7_c*xc4 )+   (exp7_d*xc3 ) +(exp7_e*xc2 )+   (exp7_f*xc )+    exp7_g  ;
    h1=(exp8_a*xc7 )+(exp8_b*xc6 )+ (exp8_c*xc5 )+   (exp8_d*xc4 ) +(exp8_e*xc3 )+   (exp8_f*xc2 )+   (exp8_g*xc )+ exp8_h  ;
    i1=(exp9_a*xc8 )+(exp9_b*xc7 )+ (exp9_c*xc6 )+   (exp9_d*xc5 ) +(exp9_e*xc4 )+   (exp9_f*xc3 )+   (exp9_g*xc2 )+  (exp9_h*xc )+ exp9_i  ;
    
    
    p=[a1; b1 ; c1; d1; e1; f1; g1; h1 ; i1    ];
   
    f = @(x)((p(1)*x.^8)+(p(2)*x.^7)+(p(3)*x.^6)+(p(4)*x.^5)+(p(5)*x.^4)+(p(6)*x.^3)+(p(7)*x.^2)+(p(8)*x.^1)+(p(9)*x.^0));
    answ=processRoots_focalonly(fminbnd(f, ycstart, ycend));
    
    if(answ>minans && answ<maxans)
        solzarray=[solzarray ;answ];
        xfarray= [xfarray ; xarray(i)];
        
    end
    
end
xarray=xarray';


end


function [sl]=processRoots_focalonly(p)

sl=p;

end
%
% this is how you generate this formula
% syms f1 f2 f3 f4 f5 f6 f7 f8 f9 foc1 foc2 xc yc real
% K=[foc1 0 xc ; 0 foc2 yc; 0 0 1]
% F=[f1 f2 f3 ; f4 f5 f6 ; f7 f8 f9 ]
% G=(K')*F*K
% GGT=G*(G')
% v=(2*trace(GGT*GGT))-((trace(GGT))*(trace(GGT)))
% vs=simplify(v)
% R = collect(vs ,foc1) to collect terms for foc1

%2*(foc1^2*(F(3,1) + F(1,1)*xc + F(2,1)*yc)^2 + foc2^2*(F(3,2) + F(1,2)*xc
%+ F(2,2)*yc)^2 + (F(3,3) + F(1,3)*xc + F(3,1)*xc + F(2,3)*yc + F(3,2)*yc +
%F(1,1)*xc^2 + F(2,2)*yc^2 + F(1,2)*xc*yc + F(2,1)*xc*yc)^2)^2 +
%4*(F(1,1)*(F(3,1) + F(1,1)*xc + F(2,1)*yc)*foc1^3 + F(1,2)*(F(3,2) +
%F(1,2)*xc + F(2,2)*yc)*foc1*foc2^2 + (F(1,3) + F(1,1)*xc +
%F(1,2)*yc)*(F(3,3) + F(1,3)*xc + F(3,1)*xc + F(2,3)*yc + F(3,2)*yc +
%F(1,1)*xc^2 + F(2,2)*yc^2 + F(1,2)*xc*yc + F(2,1)*xc*yc)*foc1)^2 +
%4*(F(2,1)*(F(3,1) + F(1,1)*xc + F(2,1)*yc)*foc1^2*foc2 + F(2,2)*(F(3,2) +
%F(1,2)*xc + F(2,2)*yc)*foc2^3 + (F(2,3) + F(2,1)*xc + F(2,2)*yc)*(F(3,3) +
%F(1,3)*xc + F(3,1)*xc + F(2,3)*yc + F(3,2)*yc + F(1,1)*xc^2 + F(2,2)*yc^2
%+ F(1,2)*xc*yc + F(2,1)*xc*yc)*foc2)^2 - (foc1^2*(F(1,3) + F(1,1)*xc + F(1,2)*yc)^2 + foc1^2*(F(3,1) + F(1,1)*xc + F(2,1)*yc)^2 + foc2^2*(F(2,3) + F(2,1)*xc + F(2,2)*yc)^2 + foc2^2*(F(3,2) + F(1,2)*xc + F(2,2)*yc)^2 + (F(3,3) + F(1,3)*xc + F(3,1)*xc + F(2,3)*yc + F(3,2)*yc + F(1,1)*xc^2 + F(2,2)*yc^2 + F(1,2)*xc*yc + F(2,1)*xc*yc)^2 + F(1,1)^2*foc1^4 + F(2,2)^2*foc2^4 + F(1,2)^2*foc1^2*foc2^2 + F(2,1)^2*foc1^2*foc2^2)^2 + 2*foc1^4*(F(1,1)^2*foc1^2 + F(1,1)^2*xc^2 + 2*F(1,1)*F(1,2)*xc*yc + 2*F(1,1)*F(1,3)*xc + F(1,2)^2*foc2^2 + F(1,2)^2*yc^2 + 2*F(1,2)*F(1,3)*yc + F(1,3)^2)^2 + 2*foc2^4*(F(2,1)^2*foc1^2 + F(2,1)^2*xc^2 + 2*F(2,1)*F(2,2)*xc*yc + 2*F(2,1)*F(2,3)*xc + F(2,2)^2*foc2^2 + F(2,2)^2*yc^2 + 2*F(2,2)*F(2,3)*yc + F(2,3)^2)^2 + 4*foc1^2*foc2^2*(F(1,3)*F(2,3) + F(1,1)*F(2,1)*xc^2 + F(1,2)*F(2,2)*yc^2 + F(1,1)*F(2,3)*xc + F(1,3)*F(2,1)*xc + F(1,2)*F(2,3)*yc + F(1,3)*F(2,2)*yc + F(1,1)*F(2,1)*foc1^2 + F(1,2)*F(2,2)*foc2^2 + F(1,1)*F(2,2)*xc*yc + F(1,2)*F(2,1)*xc*yc)^2